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Phase-transitionlike behavior is found to occur in globally coupled systems of finite number 
of elements, and its theoretical explanation is provided. The system studied is a population of 
globally pulse-coupled integrate-and-fire cells subject to small additive noise. As the population size 
is changed, the system shows a phase-transitionlike behavior. That is, there exits a well-defined 
critical system size above which the system stays in a monostable state with high-frequency activity 
while below which a new phase characterized by alternation of high- and low frequency activities 
appears. The mean field motion obeys a stochastic process with state-dependent noise, and the 
above phenomenon can be interpreted as a noise-induced transition characteristic to such processes. 
■ Coexistence of high- and low frequency activities observed in finite size systems is reported by N. 

Cohen, Y. Soen and E. Braun[Physica A249, 600 (1998)] in the experiments of cultivated heart 
^ , cells. The present report gives the first qualitative interpretation of their experimental results. 
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• Collective dynamics of coupled dynamical elements represents a central issue of nonlinear dynamics, and has served 

■ as a subject of extensive study over the last few decades. The relevant branches include chemical reactions [2], society 

of living organisms [3-5], lasers [6,7], semiconductors [8,9], neural networks [10-12] and cardiac systems [3,13,14]. In 
most theoretical studies, however, the system size is assumed infinite. While this idealization could be valid for such 
systems like spatially extended chemical reactions, there seem to be important practical cases where the finiteness 
. of the system size should explicitly be taken into account. Two existing theoretical studies on collective dynamics 
' ^ ' of populations are mentioned here in which the finiteness of the system size plays a crucial role. Firstly, Daido [15] 
. investigated the collective behavior of an inhomogeneous system of oscillators focusing on the statistics of fiuctuations 
close to the onset of collective motion. Secondly, Pikovsky et.al [16] studied coupled bistable elements and showed 
I numerically a phase-transitionlike behavior existing only for finite systems. Apart from these theoretical studies, an 
^ interesting experiment was reported recently by N. Cohen et.al [1]. They cultivated heart cells in various population 
00 ! sizes and found that the cells exhibit system-size-dependent behavior. Some more details of their reports arc the 

■ following. The heart cells extracted from ventricles of neonatural rats were cultivated, and time series of spontaneous 
spiking activities of these cells were recorded. Initially, the heart cells were made free of mutual contact by a chemical 
treatment, but after some time they spontaneously assemble to form subgroups. These cell groups were cultivated, 

'J and the spontaneous spike activity of the individual cells was observed. One interesting feature of the experimental 
results may be stated as follows. An isolated cell shows random and slow spontaneous spike activity where the 
interspike intervals (ISI) distribute sparsely with the average ISI of 5 ~ 10 seconds. The series of spikes exhibited by 
the cells belonging to relatively small groups are composed of two components characterized by a distribution with 
double peaks. The first low-frequency component is practically the same as that obtained from the spike series for 
isolated cells. Another component seems to come from a spike series which is more regular and its ISI is shorter, i.e., 
^ i about 1 second. These low- and high frequency spikes are visited alternately in time. If the cell group is sufficiently 

■ large, the individual cells exhibit only high-frequency spike activities of good periodicity. The above experimental 
results are remarkable in that the dynamics of the individual cells depends qualitatively on the population size, and 

d still await theoretical interpretation. In particular, we would like to know the origin of the transition which occurs 
as the population size is changed and whether the mechanism involved is universal beyond the particular class of 
systems of cultivated heart cells. The goal of this paper is to provide an answer to these questions. This can be 
achieved analytically by using a simple dynamical model showing coexistence of high- and low frequency activities 
under suitable conditions. Specifically, we employ globally coupled noisy integrated-and-fire model which is often used 
in the studies of neurodynamics of the brain. The dynamics of a single cell is given by 

x.,{t)^I + £,,{t), x^{t)<l, (z = l,---,iV), (1) 

where S,i{t) is white Gaussian noise with the properties {£,i{t)) — and {£.i{t)£,j{t')) = 2D ■ SijS{t — t'), and N is the 
system size which represents the principal control parameter. At the instant when Xi reaches the threshold Xi = 1, a 
spiking or firing event is assumed to occur for this cell in such a way that Xi is immediately reset to a certain value 
/. Thus, 

a:,(ii")) = l ^ a:,(<l")+d<) = /, (n = 1, • • • , iV) , (2) 
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where is the timing of the spike, the tag (n) indicates the numbering of the spiking events. When a ccU fires, 
this will immediately cause a pulsatile stimulus on all the other cells in the population. A given cell will receive a 
sum of such stimuli coming from various cells. Assuming that the effect of each stimuli decays exponentially, one may 
conveniently introduce an order parameter r{t) by 



H*) = ^E E exp{-X{t-tf'^)), (3) 

whose effect is experienced commonly over the cells. Since the primary effect of the order parameter should be to 
lower the effective threshold for each cell to fire, it would not be unreasonable to assume the dependence of the reseting 
value of / on r like 

fir) = foil -e-^n- (4) 

Note that increasing r implies increasing / with the upper limit /q. 

Before discussing collective behavior, we show how a single cell behaves under fixed r. Equations (1), (2) and (4) 
determine how the cells behave under a given value of r. To characterize statistically the sequence of spikes generated, 
we derive a density distribution function of ISI. Because ISIs arc given by the first passage time [17,18] of the stochastic 
process given by Eq.(l), their distribution can be obtained with the standard method and takes the form 

From this distribution, spike frequency w which is defined as the inverse of mean ISI, i.e., 1/{T), becomes 

"^^^ = ^ il^fo)ifo-e-^r- (6) 

One can see that larger /smaller r corresponds to larger/smaller u or higher/lower frequency activity. The above 
results are also consistent with the experimental facts for the heart cells that when the cells are isolated the spiking 
frequency is the lowest and the distribution of the corresponding ISIs is the broadest. 

Two remarks should be given on our model. Firstly, we have taken into account the effect of mean field on the 
reseting state by assuming / to depend on r. However, our choice of the specific form of /(r) is rather arbitrary. What 
is important here is that the rate of spiking lu should be an increasing fimction of r. Under this condition, other choices 
of fir) would give qualitatively the same results. Secondly, our main goal is to understand, analytically if possible, 
some general features of the collective dynamics shared commonly by finite-size systems rather than reproducing 
precisely the experimental results for the real heart cells, so that Eqs. (1) to (4), which may not be so realistic for 
heart cells, would still be useful enough as a model to work with. 

Equations (1) to (4) are calculated numerically with some values of N, from which the distribution functions for r, 
denoted as P(r, t), are obtained. Figures 1, 2 and 3 show P(r, for the cases of N = 20, 100 and 200, respectively. 
In order to make clear the difference between the cases of N = 100 and 200, the data in Figs. 2b and 3b are plotted 
in semilogarithmic scales. For sufficiently large system size, the distribution has a single sharp peak about r = 1.0 
corresponding to coherent high-frequency activity. It seems that the distribution shows no qualitative change as N 
becomes even larger. In contrast, a remarkable change occurs for smaller A''. Figures 1 and 2 show the appearance of 
a new peak at smaller r corresponding to low-frequency activity. The double peaks implies coexistence of the different 
steady states. Actually, the coexistence of high- and low frequency states is the feature observed experimentally in real 
heart cells. Our model gives a fairly well-defined critical population size associated with a transition from monostable 
state to coexistence state. 

In what follows, our analysis proceeds in two steps. We first derive an evolution equation for the mean field 
This can be achieved by assuming that the mean field r evolves much more slowly than the evolution of the individual 
Xi so that the adiabatic approximation is applicable. This assumption implies that the random variable r obeys a 
Markov process so that the equation for P(r, t) can be written in the form of a Kramers-Moyal expansion [19] [20] as, 

(-J^)"-'"'(M)P(M). P) 

where D^"'^ is defined as 
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i?(")(r, i) = 1 lim + r) - r(t)]")|,(,)=.. (8) 

To find explicit forms of D^"\ we calculate r{t + r) — r(t) from Eq.(3) for small r under the condition that r{t) = r. 
We obtain 

r-(i + r) - r{t) = (e"^- - 1)- + ^ E E e^^(*+-*^"') 

»=1 n;t<t("><t+T 
1 " 

= (e-A-_l)^+_^„,(i,t + ^), (9) 

1=1 

where ni(t, t + r) is the number of spikes of the i— th cell during a short interval t ~ f + r. To obtain the last expression 

in Eq (9), we use smallness of r and replace e~'^(*+'^~*i ^ with 1. Since r{t) is a slow variable by assumption, no 
information on the dynamics of the individual Xi is relevant except for the timing of their spiking. The sequence 
of the spikes is expressed by a Poisson process of the mean spike rate w(r), so that the distribution for ni{t,t + r), 
denoted by Pn{n), becomes 

Applying Eq.(lO) to Eq.(9) and then to Eq.(8), we obtain 

= -Ar + lim -in) = -Xr + u)(r), 

T— >0 T 

I)(2) = 1 lini _ = (11) 

27V r^O T^^ ^ " ' 2N ^ ' 

Noting that ZP^"^ ~ 0{N^~"), and neglecting the terms less than 0{N~'^), Eq.(ll) then reduces to a Fokker-Plank 
equation [21] 



dt dr 

- _^ 

dr 



i-Xr + coir))Pir,t) - -^'^Pir,t) 



(12) 



The above equation is equivalent with the Langevin equation 



r^t)=-Xr + iv{r) + ^'^.m, (13) 

where the last term represents noise which is gaussian with the properties {^{t)) = and (C(i)C(i')) = S{t — t'). From 
this evolution equation (13) and the characteristic time scale of Xi, i.e., = {I — f{r))^/D, in which Xi can diffuse 
over the interval [/(r), 1], the condition for the adiabatic approximation to hold may be expressed as 

D ^^|Z^!:±i^|. (14) 



(i-ZW) 

It should be noted that the noise strength of the stochastic process Eq.(13) depend on r. Indeed, the fact that the 
mean field of a finite size coupled system obeys a stochastic process with state dependent noise is a main conclusion 
of the present paper. 

As the second step, we study the behavior of the mean field r using Eqs.(12) and (13). The systematic part of 
Eq.(13), i.e., — Ar + uj{r), admits a single stable steady state r = Tc as shown in Fig. 4. Thus, for N = oo the 
distribution has a single delta peak aX r = Tc- The system with smaller N behaves differently, which comes from 
the state-dependence of the noise strength. It is well known that stochastic processes with such noise exhibit phase- 
transitionlike behavior called noise- induced transition [22], which means a bifurcation exhibited by the locus of the 
extrema of the distribution function. It can be shown that Eqs.(12) and (13) exhibit such a transition. The steady 
distribution for r is obtained from Eq.(12) where we assume vanishing probability current, i.e., J = 0. Thus, 
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Psteady{r) oc exp(-$(r)) 



= exp(- 



„ >r / + UJ(S) , , , . 



), (15) 



where <& represents an effective potential. Because the extrema of P are identical with the extrema of $, they can be 
found from 

or 

N= ^ . (17) 

Equation (17) gives a bifurcation diagram for the extrema as depicted in Fig. 5, where the solid lines and dotted line 
show the loci of maxima and minima of Psteady (t) , respectively; the three vertical lines indicate the particular system 
sizes corresponding to Figs. 2, 3, and 4. The location of the maxima corresponding to each vertical lines thus obtained 
analytically is in good agreement with the numerical simulation whose results are shown in Figs.l ~ 3. It is clear 
from Eq.(17) and Fig. 5 that there exists a phase transition when N is changed. There exists a well-defined critical 
population size Nc for the transition. 

In conclusion, the collective dynamics of intcgrate-and-fire cells is studied, from which the occurrence of a transition 
at a finite system size similar to cultivated heart cells is confirmed. In order to understand the origin of the transition, 
a stochastic differential equation for the mean field is first derived. The transition may be regarded as a noise-induced 
transition peculiar to systems with finite size. Such results do not seem confined to a specific model adopted, but 
could bo observable in wide variety of noisy finite-size populations. 

Finally, a few more comments should be given. It is not intended in the present paper to reproduce experimental 
results (e.g. those by Cohen) quantitatively. Our main goal was to make clear, with the aid of a relatively simple 
model, a certain qualitative feature of the collective dynamics exhibited by finite-size populations with noise. In 
more realistic models, the coupling should be local rather than global. In fact, real heart cell interact through 
electrical coupling or gap junctions which is local, and this fact is completely ignored in the present analysis. The 
first important theory on the transitions induced by finitencss of the system size was developed by Pikovsky ct al. 
[16]. Some differences between their works and ours are the following. Firstly, they assume bistability from the outset 
for the individual elements, while no such assumption is introduced in our model; bistability appears naturally as a 
result of collective dynamics there. Secondly, we succeeded in clarifying analytically the origin of the transition as a 
noise-induce transition. 

The author thanks Y. Kuramoto for fruitful discussions and critical reading of the manuscript. He also thanks A. 
S. Mikhailov for valuable discussions. 
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FIG. 1. Distribution function for the mean field r calculated numerically from Eq.(l)~(4) with N=20; other parameters are 
7 = 0.087, A = l.Q,D = 10-^/3 = 5.0 and /o = 0.92. 
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FIG. 2. Distribution function for the mean field r calculated numerically from Eq.(l)~(4) with N=100, other parameters 
are same to Fig.l, (b) is plotted with semilogarithmic scale. 
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FIG. 3. Distribution function for the mean field r calculated numerically from Eq.(l)~(4) with N=200, other parameters 
are same as in Fig.l; in (b) the data are plotted in semilogarithmic scale. 
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FIG. 4. Systematic part of Eq.(13) as a function of r. 
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FIG. 5. Bifurcation diagram of cxtrcma of the steady distribution function P{r) obtained from Eq.(17). Solid and dotted 
curves indicate maxima and minima, respectively. Vertical lines indicate the values of N chosen in Fig.l ~ 3. 
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